### Fit baseline model specification on just post-WW2 data

date()
library("dplyr")
library("purrr")
source("backend_main.r")
sessionInfo()


load("kr_analysis_dispute.rda")
load("kr_analysis_participant.rda")
load("confirm_fit_postwar.rda")


## Set to FALSE to estimate from scratch (takes much longer)
CONFIRM_MODE <- TRUE

## Identify post-WW2 disputes
id_postwar <- data_participant[[1]] %>%
  group_by(id) %>%
  summarise(min_year = min(year), .groups = "drop") %>%
  filter(min_year > 1945) %>%
  pull(id)

## Filter data down to post-WW2 disputes
filter_postwar <- function(x) filter(x, id %in% !! id_postwar)
data_dispute_postwar <- map(data_dispute, filter_postwar)
data_participant_postwar <- map(data_participant, filter_postwar)
stopifnot(all(map_int(data_dispute_postwar, nrow) == length(id_postwar)))

## Model formulas
f_dispute <- id + war + win_a + win_b_alt ~
  polity_a + majpow_a + log1p(py_alt) + s_cinc + contig |
  polity_b + majpow_b + log1p(py_alt) + s_cinc + contig |
  log(n_states_a) | log(n_states_b)
f_participant <- id + sidea ~
  log(gdp_pwt) + log1p(irst) + log1p(pec) + log1p(tpop) + log1p(upop) +
  log1p(distance) + nuclear | log1p(pct_imports) + polity2

## Fitting process
fit_postwar <- est_structwar_ac(f_dispute = f_dispute,
                                f_participant = f_participant,
                                data_dispute = data_dispute_postwar,
                                data_participant = data_participant_postwar,
                                n_halton = 1024,
                                init = if (CONFIRM_MODE) confirm_fit_postwar else NULL,
                                scale = TRUE,
                                reltol = 1e-14,
                                iterlim = 5000,
                                printLevel = 1)

## Confirm that coefficients match
if (CONFIRM_MODE) {
  for (i in seq_along(confirm_fit_postwar)) {
    stopifnot(all.equal(coef(fit_postwar[[i]]$fit), confirm_fit_postwar[[i]]))
  }
}

## Seeds generated 2023-01-17 by random.org
seeds <- c(542550, 524282, 550035, 55694, 317684,
           742063, 378961, 760389, 212451, 45493)
boot_index_postwar <- vector("list", length(fit_postwar))
for (i in seq_along(fit_postwar)) {
  set.seed(seeds[i])
  fit <- fit_postwar[[i]]
  boot_index_postwar[[i]] <- foreach (b = 1:100, .combine = "rbind") %do% {
    sample(fit$data_dispute$id,
           size = nrow(fit$data_dispute),
           replace = TRUE)
  }
}

if (!dir.exists("results"))
  dir.create("results")
save(fit_postwar, file = "results/fit_postwar.rda")
save(boot_index_postwar, file = "results/boot_index_postwar.rda")


date()
